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ABSTRACT 


We  study  the  shape  of  an  acoustic  pulse  from  a  monochromatic 
source  on  the  SOFAR  axis  as  a  function  of  range  in  the  parabolic 
approximation.  Borrowing  a  method  of  quantum  physics,  we  work 
with  the  stochastic-mechanical  interpretation  of  the  parabolic 
wave  equation:  the  motion  of  acoustic  "particles"  along  random 
rays  is  governed  by  a  stochastic  differential  equation  of 
diffusion  type.  As  recognized  by  Williamson  in  his  study  of 
electromagnetic  pulses,  the  average  pulse  shape  is  just  the 
(suitably  normalized)  graph  of  the  probability  density  function 
of  a  particle's  travel  time  along  a  random  ray.  We  discuss  an 
approach  to  finding  this  probability  density  at  long  range,  and 
we  derive  a  formula  for  the  average  travel  time  as  a  function  of 
range;  this  formula  predicts  the  well-known  transition  from  the 
characteristic  near— range  pulse  shape  to  the  characteristic  far- 
range  shape.  Knowledge  of  a  key  parameter  in  this  formula  is 
equivalent  to  knowing  the  location  of  convergence  zones.  A 
possible  test  of  this  theory  against  experimental  data  is 
outlined.  In  the  course  of  the  paper,  we  derive  new  stochastic- 
mechanical  versions  of  Ehrenfest's  theorem  and  the  virial  theorem. 


A  STOCHASTIC-MECHANICAL  APPROACH 


TO  ACOUSTIC  PULSE  PROPAGATION 
IN  THE  DEEP  OCEAN  SOUND  CHANNEL 


Introduction 

The  effects  of  inhoraogeneities  of  the  medium  on  the 
propagation  of  pulses  of  acoustic  or  electromagnetic  energy  are 
of  fundamental  importance  for  a  wide  range  of  scientific  and 
engineering  work,  including  the  design  and  operation  of  radar  and 
sonar  detection  systems,  seismic  and  ionospheric  probes,  radio 
telescopes ,  and  communication  systems.  Describing  the  spreading 
of  a  pulse  as  it  travels,  and  the  shape  of  a  pulse  as  a  function 
of  time  at  various  ranges,  are  basic  problems.  For  example, 
pulse  spreading  limits  the  ability  of  zero-forcing  equalizers  to 
reduce  intersymbol  interference  in  communication  systems  ([1], 
chapter  9).  Astonomers  have  extensively  investigated  the 
spreading  of  radio  signals  from  pulsars  ([2],  [3],  and  references 
therein) . 

In  applied  underwater  acoustics,  the  use  of  sound  pulses 
extends  from  the  earliest  efforts  to  locate  airmen  lost  at  sea  to 
the  latest  attempts  to  map  the  fine  features  of  the  ocean  floor 
with  acoustic  probes.  As  is  well  known,  the  deep  sound  channel, 
or  SOFAR  channel,  with  axis  at  the  depth  of  minimum  sound  speed 
in  the  ocean,  acts  as  an  acoustic  waveguide,  so  that  sound  from 
explosions  near  the  channel,  and  other  strong  signals  entering 
the  channel,  can  propagate  to  distances  exceeding  3000 
kilometers.  Pulse  propagation  in  the  SOFAR  channel  therefore 


Page  -  2  - 


has,  among  others,  applications  on  a  global  scalet  the  monitoring 
of  underwater  volcanic  activity  and  earthquakes  to  anticipate  the 
potentially  catastrophic  accompanying  tidal  waves.  The 
importance  of  pulse  propagation  in  these  and  other  applications 
of  underwater  sound,  including  sonar,  has  recently  stimulated 
vigorous  theoretical  and  experimental  investigation  of  the 
phenomenon.  For  accounts  of  several  recent  projects,  see  [4];  a 
more  recent  study  is  [5], 

In  the  present  paper,  we  are  concerned  with  acoustic  pulse 
propagation  in  the  SOFAR  channel;  in  particular,  we  are 
interested  in  pulse  spreading  and  pulse  shape  as  a  function  of 
range  for  a  monochromatic  source  on  the  SOFAR  axis.  One  of  the 
first  features  of  this  phenomenon  to  be  recognized  is  the  radical 
alteration  of  pulse  shape  as  the  pulse  moves  from  the  region  near 
the  source  to  ranges  far  from  the  source.  This  transformation  of 
the  shape  of  the  graph  of  the  time-varying  acoustic  intensity 
received  down  range  on  the  sound  axis  is  well-illustrated  by  the 
two  records  of  acoustic  signals  shown  below,  which  were  first 
published  in  1949  [6]  and  have  been  reproduced  in  Brekhovskikh '  s 
book  ([7] ,  p.  504): 


(see) 


FIG.  A.  Near  Range 
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FIG.  B.  Far  Range 


These  pictures  are  to  be  contrasted  with  FIG.  C  below,  a 
schematic  rendering  of  a  theoretical  pulse  shape  for 
electromagnetic  radiation  due  to  Williamson  [2]j 


FIG.  C 


In  comparing  these  figures,  note  that  any  record  after  five 
seconds  in  FIGS.  A  and  B  is  of  background  noise. 

The  differences  in  figures  A,  B,  C  are  well-understood 
qualitatively .  At  short  range  in  FIG.  A,  the  direct  ray  along 
the  sound  channel  axis,  which  carries  most  of  the  acoustic 
energy,  arrives  before  the  refracted  rays,  which  account  for  the 
(more  or  less)  exponentially  decreasing  tail.  However,  at  long 
range,  in  FIG.  B,  the  refracted  rays,  even  though  they  have 
traveled  over  longer  distances  than  the  direct  ray,  have  done  so 
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faster  sound  speeds  for  a  long  enough  time  to  overtake  the 

direct  ray.  Hence,  we  find  a  (more  or  less)  exponential  increase 
to  maximum  intensity,  followed  by  an  abrupt  decrease  to  zero 
intensity.  FIG.  C  is  different  from  both  FIGS.  A  and  B.  The 

assumptions  made  in  deriving  FIG.  C  are  that  all  rays  are 

refracted  to  some  extent,  but  that  the  speed  of  light  is  constant 
throughout  the  medium  (as  seems  reasonable  for  the  scattering  of 
light  rays  by  interstellar  dust,  for  example).  Thus,  there  is  no 
direct,  or  unrefracted,  ray;  if  there  were  one  it  would  arrive  at 
the  extreme  left  of  FIG.  C,  but  all  the  photons  travel  along 
refracted  rays  and  are  delayed  somewhat.  There  is  a  gradual 
buildup  of  intensity  to  a  maximum,  followed  by  a  gradual,  near¬ 
exponential  decrease  to  zero  intensity.  We  see,  then,  that  the 
distinctive  pulse  shapes  of  FIG.  A  and  FIG.  B,  and  the 

transformation  of  FIG.  A  into  FIG.  B  with  increasing  travel  time 
of  the  pulse  are  accounted  for  by  the  variation  of  the  speed  of 

sound  with  depth  in  the  ocean,  and  the  focusing  of  acoustic  rays 

bY  the  ocean  sound  channel  at  the  depth  of  minimum  sound  speed  on 
the  channel  axis. 

In  this  paper  we  propose  to  investigate  FIGS.  A  and  B 
quantitatively.  Since  we  are  studying  propagation  in  the  sound 
channel,  we  work  in  the  parabolic  approximation  (e.g.,  [8]).  The 

domain  of  validity  of  this  approximation  depends  on  the 

wavelength,  but  it  is  essentially  a  long-range  approximation; 
hence,  our  methods  are  more  justifiably  applied  to  FIG.  B. 

Mathematically,  the  parabolic  approximation  is  embodied  in  the 
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parabolic  wave  equation,  which  is  identical  to  the  Schrodinger 
equation  of  quantum  mechanics.  We  therefore  expect  that  ways  of 
analyzing  and  interpreting  the  Schrodinger  equation  will  be 
relevant  to  our  problem.  For  example,  Dashen  [4,  9]  has  recently 
shown  that  Feynman  path  integrals  [10,  11]  are  a  useful  tool  in 
analyzing  wave  propagation  in  random  media  in  the  parabolic 
approximation . 

The  present  author  has  recently  initiated  a  line  of  acoustic 
research  [12]  based  on  another  interpretation  of  the  Schrodinger 
equation,  Nelson's  stochastic  mechanics.  (For  early  work  in 
stochastic  mechanics,  see  [13-15].  Two  recent  review  articles, 
with  references  to  subsequent  work,  are  [16]  and  [17];  see  also 
[18].)  Just  as  Williamson  in  his  work  cited  above  on 
electromagnetic  radiation  developed  a  picture  of  a  flux  of 
photons  along  stochastic  ray  paths,  we  work  with  the  stochastic 
mechanical  picture  of  acoustic  "particles"  (or  "bundles  of 
energy")  traveling  along  random  rays  from  source  to  receiver. 
This  picture  is  especially  apt  for  analyzing  pulse  shapes 
because,  as  Williamson  recognized  and  exploited,  with  proper 
normalization  the  average  pulse  shape  is  just  the  graph  of  the 
probability  density  function  of  the  travel  time  along  the  random 
rays.  According  to  [3],  Williamson  achieved  the  first 
theoretical  prediction  of  average  shapes  of  scattered 

electromagnetic  pulses  by  determining  this  probability  density 
function . 
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A  main  goal  of  our  research  has  accordingly  been  to  derive 
the  analogous  acoustic  result,  i.e.,  to  find  the  probability 
distribution  of  the  travel  time  along  random  rays  in  the 
stochastic  mechanical  picture;  or,  more  modestly,  to  find  the 
asymptotic  distribution  as  the  range  approaches  infinity.  The 
graph  of  this  probability  density  function  (suitably  scaled,  and 
reflected  about  the  vertical  axis  to  take  into  account  the 

overtaking  of  the  direct  ray  by  the  refracted  ones)  would  then  be 
the  average  of  records  like  FIG.  B  in  the  far-field  region  (whose 
demarcation  would,  as  usual,  depend  on  the  wavelength).  We  have 
so  far  not  achieved  this  goal.  One  possible  approach  to  this 
problem  is  discussed  at  the  end  of  this  paper. 

We  do  present  in  this  paper  a  formula  (eq.  (47)  below)  for 

the  average  travel  time  along  random  rays  that  predicts  the 

transition  from  pulses  like  FIG.  A  to  those  like  FIG.  B  as  the 

range  increases.  More  precisely,  if  t  (z)  is  the  time  for  the 

direct  ray  to  travel  to  range  z,  and  E(t(z))  is  the  travel  time 

to  range  z  averaged  over  all  random  rays,  then  we  derive  a 

formula  for  E(t(z))  that  indicates  that  E(t(z))  >  t  (z)  for  small 

o 

z  and  E  ( t  ( z )  )  <  tQ(z)  for  large  z.  Our  formula  is  a  simple  one 

containing  three  parameters.  Evaluation  of  these  parameters 
would  in  principle  allow  our  formula  to  be  used  for  the  location 

of  convergence  zones  (where  E(t(z))  =  tQ(z)),  a  topic  of  some 

current  research  interest  [19];  conversely,  information  about  the 
location  of  convergence  zones  could  be  used  to  test  our  theory. 
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We  wish  to  point  out  that  the  very  definition  of  the  notion 
of  travel  time  along  a  stochastic  mechanical  ray  path  presents  a 
fundamental  problem  which  itself  is  currently  under 
investigation:  how,  if  at  all,  is  one  to  define  arclength  for 
paths  which  are  nowhere  differentiable  ([20],  p.18)  and  even  fail 
to  be  of  bounded  variation  in  any  interval  ([21],  p.395)?  A 

possible  solution  to  this  problem  has  recently  been  suggested  by 
Zambrini  and  Yasue  [22],  and  we  make  use  of  their  suggestion 
below.  It  will  be  clear  below  that  the  choice  of  definition  of 
arclength  is  crucial  for  the  work  described  in  this  paper; 
moreover,  we  believe  that  it  may  not  be  an  exaggeration  to  say 
that  the  ultimate  resolution  of  the  problem  will  play  a  large 
role  in  determining  the  future  applications  of  stochastic 
mechanics.  We  mention  in  passing  that  the  equivalent  question  of 
defining  the  angle  between  a  ray  and  the  horizontal  direction  has 
recently  been  treated  by  Palmer  [23]  in  the  framework  of  Dashen's 
path  integral  theory. 

In  our  present  work  we  choose  a  deterministic  sound  speed 
profile  which  is  a  good  approximation  to  the  "canonical"  profile 
of  Munk  ([24];  [4],  p.32)  near  the  axis  of  the  sound  channel, 
the  random  character  of  the  rays  then  results  solely  from  their 
mutual  diffractive  interference;  indeed,  it  is  suggestive  to 
think  of  the  stochastic  mechanical  picture  as  introducing  into 
geometric  acoustics  the  diffractive  effects  retained  in  the 
parabolic  approximation  [12],  such  an  interpretation  being 
analogous  to  the  semi-classical  interpretation  of  quantum 
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mechanics.  It  will  be  seen,  however,  that  our  methods  can  be 
applied  also  in  the  case  of  a  randomly  fluctuating  sound  speed 
profile;  our  methods  are  "honest"  in  the  sense  of  Keller  [25]. 

A  sketch  of  the  plan  of  the  paper  follows.  In  part  I  we 
establish  notation  and  review  some  salient  features  of  the 
parabolic  approximation.  In  part  II  we  present  the  equations 
from  stochastic  mechanics  that  we  need.  Two  of  these, 
Ehrenfest's  Theorem  and  the  virial  theorem,  are  new  stochastic 
mechanical  versions  of  well-known  classical  and  quantum 
mechanical  principles;  we  expect  they  will  consequently  be  of 
independent  interest.  In  part  III,  we  define  the  sound  speed 
profile,  arclength  of  a  random  ray  path,  and  travel  time  along 
such  a  path.  In  part  IV  we  derive  our  formula  for  average  travel 
time.  In  part  V  we  discuss  the  prediction  of  average  pulse  shape 
at  long  ranges  and  other  possible  extensions,  refinements,  and 
applications  of  the  present  work. 

I .  Parabolic  Approximation 

Ignoring  azimuthal  variations,  we  work  in  the  two- 
dimensional  coordinate  system  shown  below,  in  which  z  is  the 
range  coordinate,  increasing  to  the  right,  and  x  is  the  depth 
coordinate,  increasing  in  the  downward  direction. 
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We  choose  x=0  to  correspond  to  the  axis  of  the  deep  ocean  sound 
channel.  In  the  parabolic  approximation,  the  evolution  of  the 
acoustic  pressure  field  is  described  in  terms  of  the  propagation 

of  a  complex-valued  envelope  jf  =  ljj(x,z)  which  satisfies  the 
parabolic  wave  equation  [8]: 


i  i!  _ 
k  3  z 


+  jx  (x)i£  , 


(1) 


where  k  is  the  wavenumber,  assumed  constant,  and 

jx  (x)  =  1  -  n(x)  ,  (2) 

n(x)  being  the  depth-dependent  index  of  refraction: 

n ( x)  =  Cq/C ( x ) ,  (3) 

with  C  ( x)  equal  to  the  velocity  of  sound  at  depth  x,  and 
Cq  m  C ( o )  equal  to  the  minimum  velocity  of  sound  on  the  sound 
channel  axis.  Notice  that  with  these  conventions,  ji  is  non¬ 
negative.  Equation  (1)  is  identical  to  the  Schrodinger  equation 
for  a  particle  of  unit  mass  if  we  identify  1/k  with  Planck's 
constant  divided  by  27T,  ji(x)  with  the  potential  energy,  and  the 
range  coordinate  z  with  the  time. 


The  replacement  of  the  ordinary  wave  equation  by  (1)  rests 
on  several  assumptions  valid  for  long-range  propagation  in  the 
sound  channel,  including  the  assumption  of  gentle  refraction 

n  5  lf  (4) 

which  means  that  ji  in  (2)  is  a  small  quantity,  and  the  small- 
angle  condition 

0  <<  1  radian, 

where  6  is  the  angle  between  a  ray's  direction  at  any  point  on  it 
and  the  positive  z-axis.  Since  the  ocean  sound  channel  acts  like 
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a  vertically-thin ,  horizontally  elongated  waveguide,  (5)  is  a 
reasonable  assumption  for  rays  which  propagate  in  the  channel  to 
large  distances.  We  have  mentioned  above  the  difficulty  of 
defining  e  in  our  theory,  but  we  will  still  find  that  (5)  has 
heuristic  value. 

The  parabolic  approximation  retains  diffractive  effects,  and 
hence  is  valid  to  significantly  lower  frequencies  than  geometric 
acoustics.  Unlike  normal  mode  expansions,  the  parabolic 
approximation  can  accomodate  range- dependent  inhomogeneities:  we 
could  assume  n  =  n(x,z)  and  =  p(x,z)  above,  corresponding  to 
the  consideration  of  time-dependent  potentials  in  quantum 
mechanics.  The  stochastic  mechanical  equations  given  below  would 
have  to  be  changed  in  a  seemingly  straightforward  way.  We  do  not 
pursue  this  possibility  further  in  this  paper. 

II.  Stochastic  Mechanical  Formulas : 

In  this  section  we  present  a  brief  account  of  some  basic 
assumptions  and  results  of  stochastic  mechanics  sufficient  for 
our  needs  in  parts  III  and  IV;  we  do  not  aim  for  completeness, 
but  refer  the  interest  reader  to  fuller  accounts  in  the 
literature  [13-15,  18].  We  do  give  derivations  from  first 
principles  of  two  results  familiar  in  quantum  mechanics  but,  as 
far  as  we  know,  new  in  stochastic  mechanics:  Ehrenfest's  theorem 
(cf.  [26],  p . 216  ff.)  and  the  virial  theorem  (cf.  [27],  p.168). 
We  also  give  a  new  proof  of  Nelson's  energy  conservation  Law 
[18],  a  comparatively  recent  stochastic  mechanical  result  that  we 
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often  use  below  in  parts  III  and  IV.  Our  argument  is  actually  a 
converse  to  Nelson's  argument  in  [18],  where  he  derives  the  basic 
dynamical  equation  (16)  below  from  the  law  of  energy 

conservation;  we  derive  the  energy  conservation  law  from  (16)  and 
other  basic  stochastic  mechanical  equations. 

It  must  be  pointed  out  that,  because  the  role  of  the  time 
variable  in  the  Schrodinger  equation  is  now  being  played  by  the 
spatial  variable  z  in  (1),  the  variables  below  with  their 
customary  stochastic  mechanical  names  of  "velocity"  and  "energy" 
are  actually  dimensionless  here,  and  so  must  not  be  uncritically 
identified  with  the  actual  velocities  and  energy  of  physical 

particles;  the  introduction  of  a  natural  time  scale  would 
eliminate  this  distinction,  but  we  have  not  needed  to  take  this 
step. 

Consider  a  particle  moving  along  a  random  ray  path 
z  x(z);  here  x(z)  denotes  the  depth  of  the  particle  at  range 

This  depth  x(z)  is  a  real-valued  random  variable;  more 

specifically,  the  random  process  z  ->  x(z)  is  assumed  to  be  a 
diffusion,  so  that  x(z)  is  governed  by  an  Ito  stochastic 
differential  equation  [28-30]: 

dx(z)  =  b ( x ( z ) , z ) dz  +  /2v  dw(z),  (6) 

where  b  is  called  the  drift,  v is  a  constant  called  the  diffusion 
coefficient,  and  z  ->  w(z)  is  standard  Brownian  motion  (also 
called  the  Wiener  process).  Suppose  f(x,z)  is  a  smooth  function 
which  does  not  increase  too  rapidly  as  x  ->  +oo.  (it  is 
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difficult  to  specify  exactly  how  rapidly — see  ([14],  p.103) — 

because  the  relevant  a  priori  estimates  for  solutions  of  the 
parabolic  wave  equation  remain,  as  far  as  we  know,  unknown; 
however,  roughly  speaking,  it  suffices  that  any  low  order  partial 
derivative  of  f  multiplied  by  any  low  order  partial  derivative  of 
the  density  p  defined  below  approach  zero  so  x  approaches  +oo. 
As  a  practical  matter ,  we  can  safely  assume  such  behavior  for  any 
acoustic  quantity  for  sound  channel  propagation.)  Following 
Nelson,  we  define  the  mean  forward  derivative  Df(x(z),z)  of 
f  ( x  (  z  )  ,  z  )  by 

Df  ( x  (  z )  ,  z  )  =  lim  e|^~(-x(z+Az)  'z+Az)  -  f(x(z),z)[p  \  (?) 

Az^o  ^  z  I  z  J 

where  E{Y|Pz}  denotes  conditional  expectation  of  the  random 

variable  Y  with  respect  to  the  past  P^  of  the  x~process  at  range 

z;  i.e.,  the  best  estimate  of  Y  we  can  make  from  information 
about  x(s)  for  0 <x<z . 

The  time  reversal  of  a  diffusion  process  is  again  a  process 
of  the  same  type,  so  we  have 


dx(z)  =  b*(x(z),z)dz  +  /2v  dw(z)  , 


(8) 


for  a  suitable  backwards  drift  b* ;  and  we  define  the  mean 
backward  derivative  D*f(x(z),z)  of  f(x(z),z)  by 


D*f(x(z),z)  =  lim  e/-(x-z)  - -(-x(  2  Az  }  ' z  Az  l  |  F  \ 

az|o  i  r*j  - 

where  E{y|Fz}  denotes  conditional  expectation  of  the 

variable  with  respect  to  the  future  F  of  the 

z 


(9) 

random 


x-process  at  range 
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Now,  it  follows  from  Ito's  lemma  ([29],  p.33)  that 


Df ( x ( z ) ,  z ) 


and 


jai 

3* 


( x ( z ) ,z)  +  b(x(z) ,z)g“(x(z) ,z)  +  v~ — | (x(z) ,z)  (10) 


D*f(x(z),z)  -  ||(X(Z) ,z)  +  b* (x(z) ,z)||(x(z) ,z)  -  v^—|( x ( z ) , z ) .  (11) 

3x^ 

In  particular,  for  f(x,z)  =  x,  we  get 


and 


Dx(z)  =  b ( x ( z ) , z ) 


(12) 


D*x ( z )  =  b*(x(z),z).  (13) 

Following  Nelson  [13,  14],  we  define  the  stochastic  acceleration 
a(x(t) )  by 

a(x(t) )  =  *{DD*  +  D*D}(x(t)).  (14) 

Suppose  also  that  we  are  given  a  potential  energy  function  |i(x), 
and  so  a  force  F(x)  =  -  Then  a  short  computation  using  (10)- 

(14)  shows  that  Newton's  Second  Law  (recall  we  are  taking  m=l 
whenever  a  mass  m  of  a  particle  appears) 


F ( x ( t ) )  =  a(x(t ) ) 


is  equivalent  to 


-(v. 


2 

2 

3X 


u^) 

3X 


+  ( 


3v 
3  z 


v^) 

3x 


3x 


=  0, 


where 


u  =  (b-b*)/2 


is  called  the  osmotic  velocity  (for  historical  reasons 
[ 14 ] ,  p. 21 ) ,  and 


(15) 


(16) 


(17) 

see 


v  =  (b+b*)/2 


(18) 
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is  called  the  current  velocity.  In  equations  ( 16 ) — ( 18 ) ,  we  have 
neglected  to  indicate  that  all  functions  are  evaluated  at  x(z) 
and  z  (except  ji  =  ji(x(z)),  which  we  are  assuming  has  no  explicit 
dependence  on  z);  for  the  sake  of  legibility,  we  shall  often  make 
the  corresponding  omissions  below  when  no  misunderstanding  can 
result.  Equation  (16),  our  stochastic  version  of  Newton's  second 
law,  is  the  basic  dynamical  equation  of  stochastic  mechanics. 


There  is  also  a  basic  kinematical  equation,  which  follows 
from  the  forward  and  backward  Fokker-Planck  equations  and  the 
fact  that  (7)  and  (8)  are  time  reversals  of  one  another: 


It  also  is  a  kinematical  fact  that 


(19) 


u  =  -2-  , 

p  3x 


(20) 


where  p(x,z)  is  the  probability  density  function  of  x(z) 
Assuming  p  never  equals  zero,  (20)  can  be  written 

3R 


u  =  2  v 


3* 


(21) 


for  R  =  i  log  p.  Also  write  v  as  a  gradient:  choose  a  function 
S  so  that 


„  _  ,  2S 

V  —  zy— ■ 

3x 


(22) 


If  we  now  define  if  by 

nr  R+  i  S 

I  =  e  , 

and  choose 


(23) 


v  =  ^  ,  (24) 

then  a  simple  calculation  using  (16)  and  (19)  shows  that  f 
satisfies  the  parabolic  wave  equation  (1).  Conversely,  if  if 
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satisfies  (1),  and  R,  S,  u,  v  are  defined  by  (21)  —  (23)  ,  then  u 
and  v  satisfy  (16)  and  (19).  In  this  sense,  the  stochastic 
mechanics  embodied  in  (16)  and  (19)  is  equivalent  to  the 
parabolic  wave  equation. 

We  also  have  an  important  formula  relating  the  ordinary  z- 
derivative  of  average  quantities  to  averages  of  their  mean 
forward  and  mean  backward  derivatives  ([14],  p.98;  [15],  p.204): 

if  f ( x , z )  and  g(x,z)  are  smooth  functions,  not  increasing  too 
rapidly  as  x  ->  too,  and  if  "E"  denotes  averaging  over  all  random 
rays  z  x(z) ,  then 

^-E  (f(x(z),z)g(x(z),z)) 

=  E[ (Df (x(z) ,z) ) (g(x(z) ,z) )  +  ( f ( x( z ) , z ) ) (D* (g (x( z ) , z ) ) ) ] .  (25) 

In  particular,  with  g=l  this  equation  yields 

^E(f  (x(z)  ,z)  )  =  E(DF(x(z)  ,z)  )  ,  (26) 

and,  interchanging  f  and  g  and  putting  g=l,  we  find 

^E(f (x(z) ,z) )  =  E(D*f (x(z) ,z) ) ,  (27) 

If  we  average  (26)  and  (27),  and  then  use  (10),  (11)  and  (18),  we 

find 

d  "F 

^E(f(x(z),z))  =  E(-^(x(z),z)  +  v(x(z),z)||(x(z),z)).  (28) 

Similarly,  if  we  subtract  (27)  from  (26),  and  then  use  (10),  (11) 

and  (17),  we  obtain 

2 

E  (  U  (  X  (  Z  )  ,  Z  )f~(  x  (  Z  )  ,  Z  )  +  v  1-(  X  (  z  )  ,  z  )  )  =  0, 

9x 


(29) 
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and,  replacing  by  f,  we  have 

dX 

E(u(x(z) ,z)f (x(z) ,z)  + v  |4(x( z ) , z ) )  =  0  .  (30) 

0  x 

We  now  use  (28)- (30),  together  with  the  basic  equations  (16)  and 
(19),  to  derive  Ehrenfest's  theorem,  the  virial  theorem,  and  the 
law  of  energy  conservation.  We  shall  often  suppress  the 
arguments  x(z)  and  z  below,  as  in  equation  (16)  ff. 

Taking  the  average  of  (16)  and  using  (29)  with  f=u  yields 


E(|^  +  v  — )  +  E(^)  =  0  , 

g  2  v  g  X  c)X  r 

or,  by  (28)  with  f=v, 

£e(v)  =  E(-  j£>. 

But  (28)  with  f=x  yields 

^E(x)  =  E(v)  ,  (31) 

so  we  obtain  Ehrenf est 1 s  Theorem: 

d2  a 

— 2E ( x ( z ) )  =  E ( -  |£(x(z>)>.  (32) 

dz 

This  equation  says  that  the  average  ray  z  E(x(z))  feels  the 
average  force,  which  is  to  be  distinguished  from  -  |^(E(x(z))), 
the  force  along  the  average  ray.  These  forces  are  distinct,  and 
the  average  path  is  consequently  not  a  classical-mechanical  path, 
unless  the  function  |£  is  affine;  i.e.,  unless  the  potential  y  is 
at  most  a  quadratic  function  of  x. 


To  derive  the  virial  theorem,  begin  by  multiplying  (16)  by 
and  averaging: 

2 


-E( Vx- — ^  +  ux|£)  +  E(x||)  +  E(xv||)  +  E(x|H)  =  0, 
3x  az  ax  9x 


(*) 


X 


Now  (29)  with  f=xu  yields 
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E(ux-^i  + 
9* 


u 


+  \JX' 


A 


3x 


0, 


and  (30)  with  f=u  gives 


2vE(||)  =  -2E( u2 )  ; 

substituting  the  second  of  these  two  equations  into  the  first,  we 
find 


E(vx‘2_t  +  =  E(u2)  / 


9* 


9x 


and  using  this  in  (*),  we  get 


E(x|f  +  vx^)  =  E(u2)  -  E(x^)  . 


3  z  v~ 9x 
However,  (28)  with  f=xv  implies 

E(x^  +  vx~— )  =  ^ 


9 z  ■  -ax'  -  d^E(xv)  -  E(v  >  * 


so  comparing  the  last  two  equations,  we  have  the  virial  theorem: 


dzE(xv)  —  E(u2  +  v2)  -  E(x‘4±) 

9  x 


(33) 


rom  (28)  applied  to  f(x)  $x  ,  we  see  that  we  can  write  the 
virial  theorem  in  the  alternate  form 


^2E(ix2)  =  E(U2  +  V2)  -  Elx-jjH) 

OX 


dz 


(34) 


Finally,  we  derive  the  energy  conservation  law: 


E(^u2  +  iv2  +  ji)  =  K  , 


(35) 


when  K  is  some  constant  (the  total  energy).  We  must  show  that 


£e<}u2  +  iv2  +  y,  =  0  , 
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or ,  by  ( 28 ) ,  that 


E(if +  ^ +  - 


,  3u. 

-E  (  Urr— -  +  UV15— ) 
d  z  3x 


(&) 


We  will  prove  (&)  by  showing  that  each  side  of  this  equation 
.  3u  3v 

equals  -VE(^  ^)  .  Examining  the  average  of  equation  (16) 
multiplied  through  by  v,  we  find  that  the  left-hand  side  of  (&) 
equals 

E(vvl~f  +  uv^>  '* 

3x^  9x 

applying  (30)  to  the  second  term  transforms  this  expression  into 


vE(v^|  -  ijfvfH,, 

which  is  just  -vE(|^  |^) ,  as  claimed. 

Turning  to  the  right-hand  side  of  (&),  we  see  that  the 
average  of  (19)  multiplied  through  by  u  shows  it  to  be  equal  to 

E<vu^f  +  u2S>  •• 

3x^  3 

applying  (30)  to  the  second  term  yields 


vE(u 


32v 


3x" 


-i-(U3Z)  )  -  _,)P/  3u  3 v. 

3x  u 3x  -  vE(a3F  a?5  • 


Therefore  (&)  is  proved,  and  consequently  so  is  (35). 
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Travel  Time ,  Arclength ,  Index  of  Refraction 
If  one  tries  to  use  the  usual  expression 

t<Z>  =  »:  cu!z!>  -  e-  <:  n(x(z))do(Z! 

o 


(36) 


to  define  the  time  t(z)  required  for  a  particle  traveling  along 
the  random  ray  x(z)  to  go  from  range  o  to  range  z,  one  is 
confronted  with  the  problem  of  giving  meaning  of  the  element  of 
arclength  do-  for  the  nondif f erentiable ,  nonrectif iable  ray.  Of 
course,  we  should  have 


( do- )  2  =  (dx)2  +  (dz)2  ,  (37) 

but  the  usual  formula  dx  =  (^)dz  is  inapplicable  because  the 
•  •  dx 

derivative  doesn't  exist.  We  solve  this  problem  in  the  manner 
of  Zambrini  and  Yasue  [22],  who  have  recently  considered  this 
situation  and  proposed  the  definition 


(dx)2  =  iH^l(clz)2  , 


(38) 


where  L  is  a  certain  positive  constant.  They  choose  the  value 

L=2  by  requiring  that  (dx)2  in  (38)  have  its  ordinary  meaning  in 

t’he  limiting  case  of  non- random  differentiable  trajectories; 

however,  we  prefer  to  leave  L  as  a  parameter  in  the  theory.  In 

our  case  the  small  angle  condition  (5)  would  suggest  that  (dx)2 
2 

<<  (dz)  ,  so  that,  because  of  (35),  we  should  require 

I,2  »  K  .  (39) 


Using  formulas  provided  below,  it  would  seem  feasible  to 
estimate  the  parameters  K,  L  from  empirical  data;  condition  (39) 
would  then  provide  a  check  on  the  validity  of  our  theory. 
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It  remains  to  specify  the  index  of  refraction  n,  or 
equivalently  the  quantity  ji=l-n.  According  to  [4,  p.32],  a 

"canonical"  model  ocean  due  to  Munk  [24]  provides  a  "reasonable 
description  of  an  oceanic  sound  channel": 

LUX)  =  , 

1+6 ( ex-x-l ) 

where  6  is  a  small  dimensionless  parameter  (6=  5.7  x  10-3)  and 
the  depth  x  is  measured  in  suitable  units  (approximately  half¬ 
kilometers)  .  Since  6  is  small,  and  x  is  small  for  propagation  in 
the  sound  channel,  we  neglect  the  term  with  6  in  the  denominator, 
so  that 


j-i  (  x )  =  6(ex-x-l)  ; 

using  the  smallness  of  x  again,  we  approximate  this  function  by 
the  first  term  in  its  Maclaurin  series: 


ji  ( x)  =  ^x2  . 

Using  (37)  and  (38)  in  (36), 
.  .  z 

t(z)  =  )Q(l  -  ji(x(z)))(l  + 


we  obtain 

u2(x(z)  , z)+v2(x(z)  ,z)U-, 

2 - >  dz 


(40) 


Recalling  (35)  and  (39),  using  the  approximation  (1+w)^  =  1+iw 
for  small  w,  and  dropping  the  product  of  two  small  terms,  we 
arrive  at 


t(z)  =  fZfl  -  ^(xfz))  +  — - ^ x ( 2 ) , z ) +v2 ( x ( z ) , z  A 

o  °  2L2  2 


dz 


(41) 


Equation  (41)  serves  as  our  definition  of  the  travel  time  from 
range  0  to  range  z  along  the  random  ray  x(z).  Since  the  travel 
time  tQ(z)  for  the  direct  ray  is  just  tQ  ( z )  =  z/c0>  if  we  denote 
by  tD(z)  the  difference  tQ(z)-t(z)  between  these  travel  times. 


we 
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have  from  (41),  with  jj  given  by  (40), 

tp(z)  =  j  (  z  )  "  — ^  (u^  (x  (  z )  ,  z  )  +  v2(x(z),z))]dz  .  (42) 


2L 


IV.  Average  Travel  Time 

Let  us  take  the  average  of  equation  (42)  above;  we  get 


E(tD}  =  (f"  *o(E(!x2)  “  -4E(iu2+iv2) )dz  , 
o  n 


which,  by  (35)  and  (40)  becomes 

2 

E(tD)  =  f  (E(u)  -  -—(K-E(u) )dz  , 
o  h 

or 

E(tQ)  =  (f-  +  ( 1+1/L2  )  E(ja )  }  dz  . 

o 


(43) 


To  perform  the  integration,  we  need  a  formula  for  E(ji)  ;  this  is 
readily  obtained  from  the  virial  theorem  (34),  which  for  jj  given 
by  (40)  may  be  written 
2 

-^r-E(u)  =  6E(U2+V2)  -  26E  ( u )  . 

dz^  1  ' 

Using  (35)  again,  we  find  the  last  equation  becomes  simply 
d2 

-^■EQjO  +  46E(u)  =  26K  . 


dz' 


(44) 


we 


Supposing  the  acoustic  source  to  be  at  the  origin  of  our 
coordinate  system,  and  noting  that  ^E(jji)  =  GE(vx)  ,  by  (28), 
have  the  initial  conditions 
dE  ( u ) 


dE  ( u )  „ .  .  . 

3?  V  "  E(>>)  =  0 


z=0  . 


(45) 
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The  unique  solution  to  (44)— (45)  is 

K, 


E  ( |i ) 


•[  1  -  cos  (  2./Gz  )  ] 


(46) 


Substituting  (46)  into  (43),  performing  the  integration,  and 
recalling  tQ(z)  =  z/Co,  we  obtain 


Kt  (z) 


t  (z)-E(t(z))  =  E(tn)  = 
o  D 


o'2'  t  1 
2  ,2 


+  [1  - 

,  s  i  n  (  2  /  G  z  )  ,  ^ 

C  2/Gz  ’  Jj 

,  ( 

47) 

i  the 

introduction 

for 

the 

t  ( z ) 
o 

of  the  direct 

ray 

to 

range  z  and  the  average  travel  time  of  all  rays  to  range  z. 


Let  us  observe  some  features  of  formula  (47).  First,  note 
that  E(tD)  is  proportional  to  t  ,  in  accord  with  the  well-known 
rule  that  pulse  spreading  in  the  sound  channel  is  proportional  to 
range  ([31],  p.102).  Also  observe  that 
(t  (z)-E(t(z) ) ) 

lim  (  t  (z)  }  =  772  <  0  ' 

z->0  O  2L 

so  that  the  direct  ray  has  shorter  than  average  travel  time  and 
therefore  arrives  early,  as  in  FIG.  A  above,  for  relatively  small 
values  of  the  range  z.  On  the  other  hand,  provided 


L2  >  1  ,  (48) 

we  see  that  t  (z)-E(t(z))  >  0  for  large  z,  so  that  refracted  rays 
tend  to  arrive  before  the  direct  ray,  as  in  FIG.  B.  At  some 
intermediate  range  tQ(z)  =  E(t(z)),  and  we  have  a  convergence 
zone.  Assuming  G  is  known  from  properties  of  seawater  (as  is  the 
case,  at  least  approximately,  according  to  [4]),  equation  (47) 
implies  that  knowledge  of  the  location  of  the  convergence  zone  is 
equivalent  to  knowledge  of  the  value  of  L.  To  check  the  validity 
of  the  theory,  we  could  empirically  determine  the  location  of  the 
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convergence  zone  for  pulse  propagation  along  the  sound  channel 
axis  and  use  it  to  determine  L  from  (47)  (note  that  such  an  L 
would  automatically  satisfy  (48)).  We  could  then  use  suitably 
normalized  and  averaged  records  of  the  type  in  FIG.  B  to  form  an 
empirical  estimate  of  the  probability  density  function  of  t(z), 
and  hence  of  E(t(z)).  Knowing  6  and  L,  we  could  then  estimate  K 
from  equation  (47)  and  see  whether  (39)  holds. 

Locating  the  source  off  the  sound  channel  axis,  at  depth 

'  ^^Q^ires  initial  conditions  other  than  (45)  :  one  condition 

g  2 

is  E(^)=-  xq  at  z=0 ,  but  specifying  E(vx)  at  z=0  seems  to  require 
the  introduction  of  another  undetermined  parameter  into  the 
initial  conditions,  and  hence  into  the  (readily  obtained) 
analogue  of  equation  (47). 

Note  that  (46)  implies  that 

1\m  \  foE(j1(x(z)l)az  -  |  -  (49) 

z-»oo  *• 

which,  interpreted  as  a  quantum-mechanical  statement,  is  the  well- 
known  equipartition  of  energy  for  the  harmonic  oscillator. 
Comparing  equation  (46)  with  the  corresponding  equation  in  the 
usual  quantum-mechanical  proof  of  (49)  (cf.  [32],  p.108)  yields 

the  fact,  not  immediately  obvious,  that  the  matrix  elements  of 

the  potential  energy  between  different  and  non-adjacent 
eigenstates  vanish. 

With  y  given  by  (40),  Ehrenfest's  Theorem  (32)  and  the 
initial  condition  x(0)=0  yield 


E  ( x  (  z  )  )  =  A  sin(./ez)  , 
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where,  by  (31),  A  =  ( 1/ye)  E(  v<  0 , 0 ) )  .  Thus,  if  there  is  any 
difference  between  the  amount  of  acoustic  energy  initially 
directed  upward  and  the  amount  initially  directed  downward,  so 
that  E(v(0,0))^0,  then  the  acoustic  particles  on  the  average 
execute  simple  harmonic  motion  above  and  below  the  sound  axis. 

V.  Discussion 

We  have  presented  a  theory  of  acoustic  pulse  propagation  in 
the  deep  ocean  sound  channel  in  terms  of  the  stochastic 
mechanical  interpretation  of  the  parabolic  wave  equation.  We 
have  derived  the  formula  (47)  for  average  travel  time  along  a 
random  ray  and  shown  how  it  accounts  for  the  transformation  of 
near-range  pulse  shapes  (FIG.  A)  into  far-cange  pulse  shapes 
(FIG.  B)  .  We  have  also  described  how  the  theory  can  be  tested 
against  empirical  data. 

As  mentioned  in  the  introduction,  we  would  like  to  use  the 
theory  to  predict  the  average  long-range  pulse  shape;  or, 
equivalently,  to  find  the  asymptotic  probability  distribution  of 
tD(z)  as  z-»oo.  It  follows  from  equation  (49)  and  an  easy 
generalization  of  work  of  Darling  and  Kac  [33]  to  non-stationary 
processes  that 

1  i z 

j  )0p(x(z)  )dz  -»  |  ,  z ->oo  , 

in  probability.  If  we  could  decompose  (u^+v^)  as 
2  2 

u  +  v  =  A  +  B  , 
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nd 


where,  in  some  suitable  far-field  approximation, 

a)  A  and  B  have  no  explicit  dependence  on  z;  i.e.,  they 
depend  on  z  only  through  x(z), 

1  lZ 

b)  —  JQE(A)dz  ->  C^f  for  some  constant  C^, 

1  ( Z 

c)  L ( z )  '0E(B)dz  C2'  for  some  constant  C2  and  slowly 
varying  function  L(z)  such  that  L(z)->oo  as  z->oo  (for 
example,  L(z)=log  z,  L(z)=log(log  z),  etc.), 

then  equation  (42)  and  the  results  of  Darling  and  Kac  would  imply 

1 1  z 

that,  if  we  put  X(z)  =  -)  B(x(z) )dz ,  a  suitably  normalized 

z  o 

version  of  tD(z)  differs  from  X(z)  by  a  process  approaching  0  in 
probability  as  z->oo,  and  the  distribution  of  a  suitably 
normalized  version  of  X(z)  converges  to  a  unit  exponential 
distribution  as  z-»oo.  In  view  of  FIG.  B,  this  would  seem  an 
appropriate  and  desirable  result.  We  do  not  know  whether  such  a 
decomposition  satisfying  a)-c)  is  possible. 

To  apply  the  results  of  Darling  and  Kac,  or  the  results  of 
followers  (e.g.,  [  34  38  ] )  in  analyzing  the  asymptotic 

distribution  of  functionals  of  diffusion  processes  of  the  form 

I  f(x(z))dz  , 

condition  a)  is  essential.  These  authors  also  assume  that  the 
process  x(z)  is  temporally  homogeneous  (i.e.,  stationary). 
Because  of  the  explicit  dependence  of  the  drift  b  on  z  in 
equation  (6),  this  assumption  is  not  justified  in  our  case.  In 
attempting  to  apply  results  like  those  of  [33-38],  therefore,  one 
must  always  inquire  whether  the  authors’  arguments  hold  for  non- 
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stationary  processes.  The  same  difficulties  are  encountered  in 
trying  to  apply  other  limit  theorems  for  diffusion  processes 
(e.g.,  [39],  [40]). 

It  would  also  be  desirable  to  explore  properties  of  the 
travel  time  t(z)  of  equation  (41)  with  a  more  realistic  choice  of 
Y  than  (40);  the  formula  preceding  (40)  seems  a  natural  next 
step.  With  this  choice  of  jji  we  have  so  far  not  succeeded  in 
finding  an  explicit  formula  for  E(jj)  to  replace  (46).  We  intend 
further  study  of  this  question,  and  expect  the  search  for 

explicit  formulas  for  E(j.i)  to  involve  suitable  approximations 

which  we  cannot  presently  anticipate. 

We  have  been  avoiding  the  question  of  boundary  conditions  at 
the  sea  surface  and  bottom.  Allowing  x(z)  to  vary  only  in  a 
finite  interval,  rather  than  from  -oo  to  too,  would  involve 

consideration  of  the  boundary  behavior  of  the  diffusion  process 
near  the  endpoints  (e.g.,  [20],  chapter  4)  and  the  attendant 
analytical  complications,  and  would  seem  to  be  of  limited 
physical  relevance  for  long-range  pulse  propagation  in  the  sound 
channel.  However,  certain  advantages  might  well  accrue.  For 
example,  if  the  time  the  particle  hits  the  bottom  (or  surface) 
has  an  exponential  distribution  (reference  [41],  chapter  V,  gives 
conditions  on  stationary  processes  which  guarantee  this),  then 
taking  absorption  at  the  surface  and  bottom  into  account  by 
"killing"  particles  as  they  reach  the  boundary  introduces  the 
usual  exponential  attenuation  factor  for  acoustic  intensity.  In 

addition,  it  is  conceivable  that  the  extra  a  priori  complications 
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of  a  model  with  finite  depth  would  lead  to  amelioration  of  some 
of  the  difficulties  mentioned  earlier — for  example,  those 
associated  with  the  prediction  of  pulse  shape. 

In  the  present  paper  we  have  used  methods  which  are 
unfamiliar  in  acoustics.  However,  we  believe  that  we  have  taken 
at  least  a  small  step  toward  demonstrating  that  these  methods  can 
yield  insight  into  at  least  some  acoustic  phenomena  amenable  to 
the  parabolic  approximation.  We  conclude  by  pointing  out  another 
underwater  acoustic  phenomenon  of  current  interest  to  which  we 
believe  our  methods  are  well-suited.  Recently  Dosso  and  Chapman 
[42]  have  studied  the  leakage  of  energy  from  a  shallow  sound 
channel  in  the  Pacific  Ocean  into  the  deep  sound  channel.  In  the 
parabolic  approximation,  this  phenomenon  is  identical  to 
particles  tunneling  through  potential  barriers  in  quantum 
mechanics.  Using  stochastic  mechanical  methods,  Jona-Lasinio, 
Martinelli ,  and  Scoppola  [43,  44  ]  have  recently  gained 
significant  new  insight  into  quantum-mechanical  tunneling  between 
wells  of  multi-welled  potentials.  Following  their  approach  in 
analyzing  the  situation  studied  by  Dosso  and  Chapman  promises  to 
lead  to  another  demonstration  of  the  value  of  stochastic 
mechanical  methods  in  research  in  underwater  acoustics. 
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